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ABSTRACT 

Acoustic modes of oscillation are affected by stellar activity, however it is unclear how starspots contribute to these changes. Here 
we investigate the nonmagnetic effects of starspots on global modes with angular degree f < 2 in highly active stars, and characterize 
the spot seismic signature on synthetic light curves. We perform 3D time-domain simulations of linear acoustic waves to study their 
interaction with a model starspot. We model the spot as a 3D change in the sound speed stratification with respect to a convectively 
stable stellar background, built from solar Model S. We perform a parametric study by considering different depths and perturbation 
amplitudes. Exact numerical simulations allow the investigation of the wavefield-spot interaction beyond first order perturbation the¬ 
ory. The interaction of the axisymmetric modes with the starspot is strongly nonlinear. As mode frequency increases, the frequency 
shifts for radial modes exceed the value predicted by linear theory, while the shifts for the £ = 2,m = 0 modes are smaller than 
predicted by linear theory, with avoided-crossing-like patterns forming between the m = 0 and m = 1 mode frequencies. The non¬ 
linear behavior increases with increasing spot amplitude and/or decreasing depth. Linear theory still reproduces the correct shifts for 
nonaxisymmetric modes. In the nonlinear regime the mode eigenfunctions are not pure spherical harmonics, but rather a mixture of 
different spherical harmonics. This mode mixing, together with the frequency changes, may lead to misidentification of the modes in 
the observed acoustic power spectra. 
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1. Introduction 

1.1. Motivation: Low degree steiiar osciiiations 

Starspots are the main observed features of magnetic activity 
in stars and play a fundamental role in understanding stellar 
dynamos. They have been detected in many stars through the 
modulation of light curves (Mosser & Appourchaux 2009) and 
Doppler imaging (Strassmeier 2009), the latter suggesting po¬ 
lar and high- to mid-latitude concentrations of magnetic fields. 
Garcia et al. (2010) unambiguously detected an activity cycle in 
a solar-like star with asteroseismology for the first time, using 
data collected by the Convection, Rotation &■ planetary Transits 
(CoRoT) mission (Auvergne et al. 2009). Other excellent data 
are available from the Kepler (Borucki et al. 2010) satellite and 
even more will come from the upcoming missions Transiting Ex¬ 
oplanet Survey Satellite (TESS) to be launched in 2017 (Ricker 
et al. 2014), and PLAnetary Transits and Oscillations of stars 
(PLATO), to be launched in 2024 (Rauer et al. 2014). 

The surface distribution of stellar magnetic activity, in prin¬ 
ciple, can be measured with asteroseismology. Gizon (2002) in¬ 
vestigated the challenge of spatially resolving surface magnetic 
activity in other stars. He concluded that it would be possible 
to discriminate between a polar cap distribution and equatorial 
band distributions of activity, but only for stars with an inclina¬ 
tion angle higher than 40°, activity stronger than the Sun, and 
rotation strong enough that the individual m-components in the 
acoustic spectra could be resolved. In a later study, using an em¬ 
pirical model calibrated to the Sun, Chaplin et al. (2007) tested 
the ability to measure, with asteroseismology, the extension of 
a latitudinal activity band distribution in Sun-like stars, find¬ 


ing that the best prospects for detection are for stars with mag¬ 
netic activity concentrated at low latitudes. More recently Santos 
et al. (2012) studied the indirect (nonmagnetic) effects on radial 
acoustic oscillations induced by the changes in the stellar stratifi¬ 
cation due to starspots. In the case of the Sun they found that the 
frequency changes are too small (by two orders of magnitude) 
to explain the observed shifts. They concluded that the indirect 
effects of starspots on the stellar stratification cannot be responsi¬ 
ble for the observed changes in the acoustic oscillations, at least 
for a star with a solar-like level of activity. 

1.2. Our approach: Time-domain numericai simuiations of 
waves and starspot 

In this work we extend the investigation to starspots in stars with 
a higher level of activity: those stars likely have starspots with 
larger surface coverage, thus implying considerable changes in 
the internal stratification (e.g., in the sound speed). We also ex¬ 
plore the possibility of identifying starspots through asteroseis- 
mic observations of highly spotted stars. We address two main 
questions: what are the changes in the mode frequencies and am¬ 
plitudes in such stars? And, what is the starspots seismic signa¬ 
ture in synthetic light curves? 

We simplify the problem by focusing on the interaction be¬ 
tween the acoustic wavefield and a single starspot with a polar 
geometry (such a configuration is compatible with Doppler ob¬ 
servations of a young solar analogue with strong activity, Mars- 
den et al. 2005). Eor the starspot we model only the indirect 
changes induced in the sound speed in the stellar interior. These 
changes in principle are not small, therefore, unlike in Santos 
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et al. (2012), a perturbative approach may not be appropriate, 
and we need to use direct numerical simulations accounting for 
the full 3D structure of the star. Direct 3D numerical simulations 
also provide synthetic observations, and therefore are a power¬ 
ful tool to characterize how the observable quantities change de¬ 
pending on the 3D stellar background. 

The remainder of this work is organized as follows: in the 
next section we introduce the GLobal Acoustic Spherical Sim¬ 
ulator (GLASS) code, which extends the code developed by 
Hanasoge et al. (2006) to include treatment of the center of 
the star. GLASS simulates linear acoustic waves propagating 
through the full 3D stellar interior, in the time-domain. In Sects. 
3 and 4 we describe the spot model employed, analyze the ef¬ 
fect of the spot on acoustic modes, and discuss the changes to 
the eigenmodes in the nonlinear regime. Section 5 highlights the 
signature of the starspot in the power spectra of the modeled light 
curves. Finally in Sect. 6 we summarize our findings and discuss 
the observational consequences. 


2. Numerical method 

2.1. Time-domain pseudo-spectrai simuiations in sphericai 
geometry 

The GLASS code solves the 3D linearized hydrodynamic equa¬ 
tions in a spherical domain containing the full star, from the 
center up to the stellar surface. The use of the linear approxi¬ 
mation is justified since acoustic wave perturbations in the Sun 
and solar-like stars have much smaller amplitudes compared to 
the stellar background quantities (e.g., velocity perturbations at 
the surface are < 20 cm/s, four orders of magnitude smaller 
than the local sound speed, see Libbrecht 1988). For the stel¬ 
lar model we considered a spherically symmetric static equilib¬ 
rium described by Model S (Christensen-Dalsgaard et al. 1996) 
stabilized against convection (Papini et al. 2014) and including 
the photosphere up to R = 1.0007 Rq, Rq being the solar ra¬ 
dius: that is our quiet Sun (QS) background model. We then 
added the starspot model to the background. We stress here that 
the spot model must not reintroduce convective instabilities and 
must fulfill the condition of hydrostatic equilibrium. We also ne¬ 
glected the perturbation to the gravitational potential (the so- 
called Cowling approximation. Cowling 1941), which reduces 
the order of the equations by 2. The use of the Cowling ap¬ 
proximation is not, in general, appropriate when considering a 
full stellar model including the center. However, since here we 
are only concerned with the changes in the modes of oscillation 
caused by perturbations to the stellar model in the near-surface 
layers, we expect the Cowling approximation to be reasonable. 
Finally we assumed that the waves propagate adiabatically inside 
the star. With these assumptions, the main equations are 


= V, (1) 

5,v = - — - —goCr - F(r)v, (2) 

Po Po 

where f (r, t) is the vector displacement of the wavefield, v(r, t) is 
the Eulerian velocity perturbation, go{r)er and po(r) are the un¬ 
signed (radially directed) gravitational acceleration and density 
of the stellar background, and r and Cr are the position vector and 
the unit vector in the radial direction. Here F(r) is a sponge-like 
damping term that prevents spurious waves reflection at the up¬ 
per boundary (Hanasoge et al. 2006). The symbols p(r, f) and 
p(r, f) are the Eulerian perturbations to density and pressure. 


both functions of ^ by means of the linearized continuity and 
adiabatic equations 

p = -V ■ (po^), (3) 

p = -^-Vpo-po(4 + Ac2)V-f, (4) 

where Ac^ is the change in the squared sound speed induced 
by the starspot, po and c^ = yoPolpo are respectively the pres¬ 
sure and the square of the adiabatic sound speed of the stellar 
background, and yo is the first adiabatic exponent. The wave- 
field is excited by setting an initial condition for the displace¬ 
ment then the code performs the temporal integration for ^ 
and V, while the other relevant quantities are computed at each 
timestep. 

GLASS employs a pseudo-spectral scheme, which uses 
spherical harmonic (SH) decomposition on spherical surfaces to 
compute the horizontal derivatives, and a sixth-order tridiago¬ 
nal compact scheme (Lele 1992) for radial derivatives. In SH 
space the spectral components are identified by the angular de¬ 
gree { and the azimuthal order m. Temporal integration is per¬ 
formed by means of a five-stage second-order low dissipation 
and dispersion Runge-Kutta (LDDRK) scheme (Hu et al. 1996). 
Eor a more detailed description of the code see Hanasoge et al. 
(2006). The grid size is chosen according to the desired reso¬ 
lution: the maximum value of f in the SH transform sets 
the minimum number of latitudinal and longitudinal grid points 
(nion = 2niat > 3f’max), while in the radial direction we adopted 
a grid with a step size constant in acoustic depth, small enough 
to resolve the background model and the shortest wavelength 
among the modes of interest. Upper boundary conditions are set 
by imposing a vanishing Lagrangian perturbation of the pressure 
(which translates into the condition V ■ f = 0 at r = R), while 
at the center of the star we prescribe regularity conditions, as we 
will now explain. 


2.2. Extending the simuiation to the center of the star 

Particular care must be taken when considering the center of the 
star. This point is a geometrical singularity in spherical coordi¬ 
nates, therefore, when calculating the radial derivatives for 
and p, boundary conditions at the center are imposed by consid¬ 
ering the asymptotic behavior of radial displacement and pres¬ 
sure in SH space. The spectral components and pr^ behave 
like and r^ respectively when r —> 0 (see, e.g., Unno et al. 
(1979)). We enforce these conditions at the center by requiring 
that the radial derivatives of and * vanish at the 

center. We then obtain the original derivatives dr^r.tm and drPtm 
by means of algebraic formulas. This procedure also ensures that 
the numerical accuracy of the compact scheme is preserved near 
the center. 

The horizontal spatial resolution increases with depth, ow¬ 
ing to the clustering of the grid points, while the radial resolu¬ 
tion remains roughly constant. As a consequence, a numerical 
instability appears in those spectral components for which the 
horizontal wavenumber kh.t - l)/r exceeds the Nyquist 

wavenumber kr.mux - njAr in the radial direction, Ar being the 
radial grid step. We remove this instability by employing a spec¬ 
tral filter that, at each radial grid point rp sets to zero all the 
spectral components satisfying kh.tirj) > kr.max- For a grid with 
constant Ar this translates to the condition 


j < it = 


sj£({ + 1 ) 
n 


(-H offset). 


( 5 ) 
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where j identifies the y-th grid point from the center. A safety off¬ 
set parameter is also implemented (we found that 2 is the small¬ 
est value that removes the instability). Finally, a fourth-order 
tridiagonal compact low-pass filter (Lele 1992) is applied in the 
radial direction, to avoid spectral blocking (Hanasoge & Duvall 
2007) and smooth the discontinuity in r caused by the spectral 
filter. The physical solutions are not affected by the filtering; ra¬ 
dial modes (the only modes that propagate through the center 
of the star) are not filtered, and all the other modes are already 
evanescent at the points rj < rj^ (the lower turning point for an 
{ - I mode at 5 mHz is 0.03 Rq, i.e., the tenth gridpoint in the 
grid we employ, while y)>=i = 2). However the combined action 
of the filter with the time evolution scheme introduces damping 
in the wavefield, with an exponential dependence on frequency. 
This results in a lifetime of ~ 50 days for waves with frequency 
at 2 mHz and down to ~ 3 days at 5 mHz for the simulation setup 
we used, which employs a radial grid with = 800 grid points 
and a time step of 1 s. 

2.3. Validation: Comparison with ADIPLS normal modes 

We validated GLASS through independent calculations of the 
theoretical modes of oscillations for our QS model. These modes 
are uniquely identified by three quantum numbers: the radial or¬ 
der n, the angular degree {, and the azimuthal order m. With 
the ADIPLS software package (Christensen-Dalsgaard 2008), 
we computed the eigenfrequencies wo.nr (which are degenerate 
in m in the QS model) and eigenfunctions of all the acous¬ 
tic modes in the range 3 < n < 37 and up to ^ = 2. We then 
used the eigenfunctions given by ADIPLS to excite the initial 
wavefield displacement in GLASS. Test simulations were run in 
a temporal window of 10 days (solar time), using different ini¬ 
tial conditions (by exciting either one or a few eigenmodes at 
the same time) and grid sizes. We also performed one simulation 
with all the modes excited and in a temporal window of 70 days, 
to use as reference for our study. 

In all the tests performed, the extracted eigenfrequencies 
and eigenfunctions were compared with ADIPLS solutions and 
showed good agreement: the difference between ADIPLS and 
GLASS for the eigenfrequencies was below 0.16 pHz (i.e., the 
frequency resolution of the simulations) and the maximum dif¬ 
ference for the eigenfunctions was ~ 0.1%, after the damping 
was accounted for. 

3. Frequency shifts: Nonlinear dependence on 
perturbation amplitude 

3.1. Spot model: Perturbation in sound speed 

We modeled the changes to the stratification caused by a starspot 
as a positive change Ac^ in the squared sound speed, while the 
density and pressure were unchanged. This guaranteed the com¬ 
patibility with the stabilization method used for the background 
model (for which a decrease in the sound speed would rein¬ 
troduce convective instabilities). We defined the change in the 
sound speed as 

Ac^(r, e) = ecl(r)f(r; (6) 

with positive amplitude e, a radial profile 

fir, rQ = ^ [cos (|r - rj/cr) -H 1] exp [-(r - r^f l{2cp-)\ 

for |r- rd < Tier and zero otherwise, and with a latitudinal profile 

giO) = ^ [cos {kO) + 1] 


for k 6 < n and zero otherwise, where 6 is the colatitude. The 
functions / and g have continuous derivatives everywhere and 
define a spot located at the north pole, with a Gaussian profile in 
radius multiplied by a raised cosine and a raised cosine profile in 
latitude. For the study the vertical and horizontal size of the spot 
were fixed by setting cr - 0.01 R© and k - 2.4, and we varied 
the depth Rq - and the amplitude e. We note that the choice 
of the coordinate system here is completely arbitrary, since we 
are studying the perturbation to a spherically symmetric (e.g., 
nonrotating) background, therefore the results of this work can 
be translated to a spot located at any point at the surface, via a 
rotation of the coordinate system. 

3.2. Linear theory 

We first discuss the effect of small perturbations, in order to 
determine for which parameter range the frequency shifts falls 
in the linear regime. We consider the normal mode solutions 

of our reference QS model, solving 

the wave equation 

POWo,«r^O = -V(CoPoV ■ ^o) - V(^o ' Vpo) + goBrV ■ (po^o) = ^^0 

(7) 

that can be derived from Eqs. (1-4). Since the reference back¬ 
ground is spherically symmetric, the solutions are degenerate in 
m. We consider now the perturbation ATT to the wave operator 
"H due to a generic change in the sound speed (while retaining 
pressure and density constant), given by 

A'7/^o = -V(Ac%V-^o); (8) 

the expression for the linear frequency shift (see, e.g., Aerts et al. 
2010) is 

^ f^Ac^(V-ro)(^-^o)podV 
OJOM 2coI^JJ^o\MV ’ 

where integration is performed over the entire volume of the star. 
Both the surface integral and the perturbation to the mode inertia 
do not contribute to the change in Aojnim in this case: the sur¬ 
face integral is zero because of the surface boundary conditions 
imposed for the pressure, and there is no contribution from the 
perturbation to the mode inertia, since the density of the back¬ 
ground is unchanged. Finally, by exploiting the separable form 
of the polar spot profile (6), we may write the frequency shift 
in terms of a product of separate integrals for the horizontal and 
radial coordinates. It is worth noting that, unlike the case of split¬ 
ting induced by rotation, the frequency shift depends only on \m\ 
and therefore modes with the same \m\ are still degenerate. This 
holds true also in the nonlinear case. 

3.3. Numerical simulations 

3.3.1. Initial conditions: (5-function source 

We performed our study in the parameter range 0.01 < e < 1 
and 0.97 R© < < 1 R©- The simulations run for 70 days 

(solar time) to reach the desired accuracy of ~ 0.16 pHz in the 
frequency domain. Wavefield displacement and velocity records 
were taken with a cadence of 60 seconds (solar time), mimick¬ 
ing the usual cadence of helio- and asteroseismic observations 
and in order to have a Nyquist frequency of ~ 8.3 mHz (above 
the maximum acoustic cutoff frequency of ~ 5.3 mHz). Starting 
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mode frequency (pHz) 


Fig. 1. Relative frequency shifts Av/v vs. mode frequency v for modes 
with £ = 0 (upper panel), { = 1 (middle panel) and £ = 2 (bottom 
panel), for a polar spot model relative to QS model with = 0.98 Rq 
and e = 0.4. Dashed lines represent linear frequency shifts given by Eq. 
(9). Fitted frequencies from numerical simulations (asterisks for m = 
0, triangles for \m\ = 1, and crosses for |m| = 2) show the nonlinear 
behavior of the shifts for m = 0 modes (asterisks). We note that m = 0 
and \m\ = 1 mode frequencies for £ = 2 cross at v - 3060 /rHz. 


from the ADIPLS eigenfunctions, we set the initial conditions 
for the waveheld displacement as 

37 2 f 

^(r,f = 0) = ;^2 (10) 

n-3 ^’=0 m=-i 

which excites all the modes from 3 < n < 37 and up to = 2. 
The initial velocity v is set to zero. We note here that all the 
modes were excited with the same phase. 

3.3.2. Nonlinear frequency shifts 

The simulated frequencies have been extracted by taking the 
SH transform coefficient ptmira, f) of the waveheld pressure 
p{r(), 6,4>, t) at each timestep and at ro = Rq + 200 km above 
the surface. This was followed by a Fourier transform in time to 
obtain the held p(m{ro, at) and then the power spectrum P{,„{a)) - 
\p(miro, w)p. Finally we divided each (£, m) spectrum in chunks 
with size of 80 pHz and centered on the peak closest to the origi¬ 
nal ADIPLS mode frequency. A least-squares Lorentzian ht was 
applied to extract frequency, amplitude, and half width at half 
maximum (HWHM) of each mode. 

Figures 1 and 2 show selected results obtained for the fre¬ 
quency shifts Av = Aa>l2n induced by the spot with respect 
to QS in the case of a polar spot located at = 0.98 Rq (in 
the following we will always show plots related to simulations 



amplitude e 


Fig. 2. Mode frequency shifts Av vs. spot amplitude e for modes 
with £ - 2, n = 12 and |m| = 0,1,2 in the case of a polar spot 
model with = 0.98. The QS eigenfrequency for these modes is 
Vo = 1970.50 ± 0.16 yuHz. Dashed lines represent linear frequency shifts 
given by Eq. (9), red asterisks are frequencies shifts from numerical 
simulations. Solid red lines represent parabolic fits to the shifts. The 
horizontal dotted line shows the crossing of the frequencies between 
the m = 0 and |m| = 1 modes at e 0.92. 


performed for this value of r^). Figure 1 shows the relative fre¬ 
quency shifts Av/v as a function of mode frequency v = coIln 
for all the modes we excited, extracted from a simulation with 
a spot of amplitude e — 0.4. The linear frequency shifts calcu¬ 
lated from Eq. (9) reproduce the behavior of nonaxisymmetric 
modes (m + 0). However for axisymmetric modes the frequency 
shifts deviate from the theoretical linear shifts, but with a dif¬ 
ferent behavior depending on t, the most interesting being the 
case £ — 2 (bottom panel). Here the relative frequency shifts 
for m - 0 modes are smaller than predicted, and even decrease 
with increasing mode frequency, while Av/v remains roughly 
constant for the \m\ - 1 and 2 modes. As a consequence, at a 
frequency of ~ 3000 pHz the Av/v for m - 0 modes cross 
the relative frequency shifts of the m - \ modes. Above this 
crossing frequency the m - I modes have the largest frequency 
shifts. Results from the other simulations performed shown that 
this crossing frequency decreases either when the amplitude e 
increases (as seen, e.g., in Fig. 2) or when is moved toward 
the surface, the latter indicating a stronger nonlinear response 
of the system to smaller surface changes than to bigger changes 
buried deeper in the convection zone. 

Figure 2 shows the simulated frequency shifts Av against 
spot amplitudes e and for all the modes with £ - 2,n — 12. 
These modes in the reference QS model are degenerate with re¬ 
spect to the azimuthal order m, with a frequency vq = 1970.50 + 
0.16 yuHz. Again linear theory successfully reproduces the fre¬ 
quency shift for nonaxisymmetric modes, shifts for \m\ - 1 how¬ 
ever start to deviate from linear behavior at e ^ 0.8. The fre¬ 
quency shifts for m = 0 modes on the contrary are nonlinear 
already at e values of ~ 0.2, with Av values smaller than pre¬ 
dicted by linear theory. A parabolic fit is able to model the shift, 
thus indicating that in this case a second order perturbative cor¬ 
rection could recover the actual frequencies. The figure shows 
also the crossing frequency (horizontal dotted line) matching the 
mode frequency at e ^ 0.92, with a value of about 1985 pVfz 
(corresponding to a frequency shift Av ^ 14.5 /rHz), decreased 
by ~ 1000 juHz with respect to the case with e — 0.4, as already 
noted above. 
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Fig. 3. Normalized cuts at dilferent colatitudes 9 of the radial displace¬ 
ment eigenfunction {r^6) at the center of the spot (0 = 0, solid 
black line), at the edge (9 = 60°, solid blue line), and outside the spot 
(dashed lines), of the £ = 2, m = 0, n = 12 mode (top) and I = Q, 
m = 0, n = 16 mode (bottom), for a polar spot with = 0.98 R© and 
e - 0.4. The solid red line is the corresponding from ADIPLS. The 
solid black line in the bottom panel beneath the eigenfunctions shows 
the radial profile / of the spot. The blue and black dashed lines in the 
top panel completely overlap. 

4. Perturbations to the eigenfunctions 

The polar spot introduces a nonspherically-symmetric perturba¬ 
tion to the background model. As a consequence the eigenfunc¬ 
tions depart from their original horizontal dependence and get 
mixed with other spherical harmonics. The extent of the mixing 
depends on the amplitude of the perturbation. Because of the ax- 
isymmetric profile of the spot there is no mixing between spher¬ 
ical harmonics with different m. We investigate the effect of the 
mixing in the eigenfunctions by writing the radial displacement 
eigenfunction of a given mode for the model with the spot as 


for = 2,m - 0,n - 12 (top panel) and for a radial mode with 
n = 16 (bottom panel), for a spot with e - 0.4 and = 0.98 R©- 
Each cut has been renormalized to its maximum. In the case of 
the quadrupolar mode the phase shift of inside the spot with 
respect to the QS solution increases from the center, and reaches 
its maximum at the edge of the spot (i.e., for 9 - 60°). The shift 
then decreases, and smoothly matches the QS eigenfunction at 
0 - 90°. The phase shift has a different behavior in the case of the 
radial mode with n = 16, it is maximum at the center of the spot 
and then approaches zero at the antipodes {0 - 180°). This is in 
agreement with what was found by Santos et al. (2012), except 
that in our case the radial profile of the eigenfunction does not 
match the QS eigenfunction at the spot edge. 

5. Synthetic power spectra for a polar spot 

5.1. Synthetic light curves 

Mode mixing affects observed light curves, as was already 
pointed out by Dziembowski & Goode (1996) and Cunha & 
Gough (2000) in the case of RoAp stars, since it changes the ex¬ 
pected mode visibilities. The production of a realistic synthetic 
light curve is a nontrivial task, which requires modeling all the 
contributions to the emergent intensity at the photosphere in the 
observed wavelength range. Here we only model the contribu¬ 
tion of the photospheric pressure perturbations to the intensity 
fluctuations I induced by the acoustic wavefield. This in princi¬ 
ple requires an explicit relation between the mode displacement 
and 7(0', (p', t) at the stellar surface, which is rather complicated 
in general (see Toutain & Gouttebroze 1993). 

For the sake of simplicity we assume that I{0, (p, t) is pro¬ 
portional to the Eulerian pressure perturbation p(ro, 9, cp, t) mea¬ 
sured at vq - Rq + 200 km. We then express I in the frequency 
domain as 

^max r 

1(9, (p, ai) oc p(ro, 9, cp, to) = EE pem(ro,tdi)Y"/(9,(p), (13) 

£=0 m=-C 


f 

t'max 

“(r, t) = Yj af'"(r)Y^(9, p) = “(r, 

f=0 


( 11 ) 


with the spherical harmonics T™ given by 


Yf(9,<p) 


(2e + 1) (€ - m)\ 
A-n ({ + m)\ 


Pjicos 6»)e™‘^, 


where ^max gives the spectral resolution of the SH transform and 
where we labeled with 7 the new mixed eigenmode originally 
represented by a pure Tf" in the QS model. Here f"^'”(r, 9) is the 

meridional profile of expressed through a truncated series 
of Legendre polynomials (i.e., the latitudinal components of the 
Yf) 


C\r,9) = 2 
^=0 




(27h- l)(7-m)! 
Ttt (t + m)\ 


7’™(cos 9). 


( 12 ) 


The amplitude coefficients a'f"‘(r) with { A £ give the degree of 
mixing with T"' of the n7m mode at the radial position r. 

Figure 3 shows radial cuts taken at different colatitudes 9 
(i.e., at different angular distances from the center of the spot) of 
the scaled radial displacement eigenfunction (r^p)'^^^"^'"(r, 0), 


where p(„, are the SH coefficients of the pressure, obtained as de¬ 
scribed in Section 3.3 and containing all the contributions to the 
wavefield in the spectral component (7, m), including the mixing 
from other £ A P modes. We can explicitly quantify the mixing 
of a single mode by writing (see also Eq. (11)) 

f 

‘'max 

p(ro, 9, (p, aj„(J = Y P(m(ro, (0, p), (14) 

e=o 

where is the mode frequency in the presence of the spot. 
The values of pf„,(ro, w„f„,) with I + t give the degree of the 
mixing. 

Figure 4 shows the coefficients \ptm(ra,co„i„y\ of expansion 
(14) for the n - 12, £ - 2, m - 0 mode for three different values 
of e. Mixing in the e - 0.01 case is almost absent. For e - 0.4 
the mixing contribution from the dipole component is negligi¬ 
ble (one order of magnitude smaller), while that from the radial 
component is 30% of the leading mode. The mixing increases 
significantly for e = 1, reaching 15% for 7=1 and 72% for 
7 = 0, thus resulting in extreme distortion of the eigenfunctions 
for that mode. 

5.2. Example power spectrum for£ - 2 

As an example we consider now four synthetic power spectra 
Pia>) of a star with an inclination angle i of the polar axis with re- 
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0 1 2 
angular degree t 


Fig. 4. Amplitude coefficients |pfm(ro, of equation (14) showing 
the mixing of the f = 2,n = 12,m = 0 mode for three simulations 
with spot amplitude e = 0.01 (red squares), 0.4 (blue asterisks), and 1.0 
(black triangles). 


We note that because of the initial conditions we set (i.e., the 
choice made in Eq. (10) of using the same phase in exciting all 
of the modes), the \m\ > 0 peaks in these spectra are twice as 
high as the peaks in a spectrum of acoustic oscillations resulting 
from stochastic excitation (that is the case for Sun-like stars). 
The e = 0.01 case falls in the linear regime: the visibility of the 
modes is the same as that of pure spherical harmonics (observed 
amplitude ratios correspond to Fig. 2 of Gizon & Solanki (2003), 
once the residual degeneracy in the m + 0 modes and the choice 
of the phase in the initial conditions (10) are taken into account) 
with very little mixing, as already shown in Fig. 4. 

The same holds true in the case of e = 0.4 and 1 for nonax- 
isymmetric modes (the differences in shape and height between 
the peaks in the plots occur because the bin size is comparable 
to the HWHM of the peaks for these modes). The m - 0 mode, 
on the contrary, departs from the linear regime both in frequency 
shift (which for e = 1 becomes even smaller than that one for 
\m\ - 1) and observed amplitude, as a consequence of the mix¬ 
ing with other spherical harmonics. 
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Fig. 5. Acoustic power vs. Av for f = 2, n = 12 and m = 0 (red) 
, ±1 (violet), and ±2 (blue) at an inclination angle i - 80°, for QS (e = 
0, upper panel) and for a spot located at a depth of 0.98 R© and with 
e = 0.01 (second upper panel), 0.4 (third upper panel), and 1.0 (bottom 
panel). The black dashed lines represent the spectra resulting from the 
sum of all the m components. Because of the initial conditions in Eq. 
(10), the observed acoustic power in the \m\ > 0 modes is twice the 
power that would be observed in the case of stochastic excitation. 


spect to the line of sight of 80°, for QS and three different ampli¬ 
tudes with e values of 0.01, 0.4, and 1, respectively. To create the 
spectra we followed the procedure outlined by Gizon & Solanki 
(2003): starting from the intensity fluctuations approximated by 
Eq. (13), we calculated the disk integrated intensity /(w) in the 
frequency domain, accounting for projection and limb-darkening 
effects (Pierce 2000). The resulting power spectra P{uj) - |/(w)P 
vs. Av are shown in Fig. 5 for the { - 2,n - \2 multiplet. 


6. Conclusion 

Using 3D linear numerical simulations, we investigated the 
changes in global acoustic modes with I <2 induced by a local¬ 
ized sound-speed perturbation with relative amplitude e, mim¬ 
icking the changes caused by a starspot with a polar cap config¬ 
uration. 

The interaction of the wavefield with a polar spot strongly 
affects the axisymmetric modes, which show a nonlinear behav¬ 
ior increasing with e, for e > 0.2. The frequency shifts for ra¬ 
dial modes exceed the shifts predicted by linear theory, while 
the shifts for the ( - 2,m - 0 modes are smaller than predicted 
by linear theory and cross the \m\ - 1 modes at a frequency that 
decreases with increasing e. For modes with m 4 0 linear theory 
successfully predicts the correct frequencies. 

The nonlinear changes, with e, in the mode frequencies and 
mode mixing (resulting from the distortion of eigenfunctions) 
will play a role in the correct identification of the modes. Strong 
mode mixing may also cause { > 2 modes to become visible in 
the observed spectrum (see, e.g., Dziembowski & Goode 1996). 
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